AD-A208  183 


rn^  r 


NAVAL  POSTGRADUATE  SCHOOL 

Monterey,  California 


DT!C 

ELECTE 
MAY  3  0 1989 


% 


0 


ACTIVE  DAMPING  OF  VIBIL^TION  IN  LARGE 

SPACE 

STRUCTURES  USING  A  KARHUNEN-LOEVE  RE¬ 
DUCED  ORDER  MODEL 

by 

Terence  M.  Grogan 
March  1989 


Thesis  Advisor  JefTB.  Burl 


Approved  for  public  release;  distribution  is  unlimited. 


’.a  Reror:  Scca-iiv  C:assa':ca;;.'r'.  L  nclas^if^ed 


RtPORI  DOCTMEMAllON  PAGE 


lb  Resinciive  Markings 


'n  Rcrrr:  N 


'a  Name  of  Periori;!;:’!  Organ;?a'.ion 
Na\al  Postgraduate  School 


■:c  Address  iv.r.  sraa’  u'lJ  ZIP  coJc  j 

Monterey.  C.A 


sa  Naa'.e  of  Funumg  Sponsoring  Organization 


6b  Office  Symbol 
r  ;r'  applkabie  i  33 


8b  Office  Symbol 
I  if  appUcabU  > 


'c  Address  <  cir.-. 


r.j  ZIP  tOJV  ' 


3  Distribution  Availability  of  Report 

Approved  for  public  release;  distribution  is  unlimited. 


5  Monitoring  Organization  Report  Numberlsi 


7a  .Name  ol  .Monitoring  Organization 
Naval  Postgraduate  School 


/b  .Address  icity.  state.  anJ  ZIP  code! 

Monterey .  C.A  93943-5000 


9  Procurement  Instrument  identification  Number 


10  Source  of  \'*rr.vNs'r< 


Program  Element  No  Project  Nc  M  ask  No  Work  L  nit  Accession  \o 


;;  rule  A-rv  se  ar;.-.  dassiCatk^o  ACl  IVE  DAMPING  OP  VIBRATION  IN  LARGE  SPACE  STRI  C  EURES  USING  A 
K.\RHUNEN-I  OEVE  REDUCED  ORDER  MODEL 


’orsc-al  Author  ?'  Terence  M.  Grocan 


;  ,ta  type  o:  Report 
Ma'ter  s  Thesis 


ih  Supplem.entary  Notation  1  he  views  e.xpressed  in  this  thesis  are  those  of  the  author  and  do  not  reflect  the  oflicial  policy  or  po- 
'ition  of  the  Deptulment  of  Defense  or  the  U.S.  Government. 


Co -at:  C  j'.'  I  IS  Subject  1  erms  i  continue  on  reverse  if  necessary  and  identity  ty  tlock  nwr-.i'er ' 

Karhunen-Loeve  Expansion.  Reduced  Order  Model.  Large  Space  Structure,  vibration 
damping. 


i.tc  ■[ ime  Covered 
From  To 


Id  Date  of  Report  lyear,  month,  day 

March  19S9 


1 5  Paze  Count 
72 


Ih  -yhstract  i  ccr.tirrut  '•'!  reverse  if  ncu'rsary  ar.d  identify  by  block  number  i 

Large  sp.ice  structures  are  dilTieult  to  control  because  of  the  higli  order  of  their  mathematical  models.  The  higJi  order 
mathematical  model  makes  the  use  of  a  reduced  order  model  to  control  the  structure  desirable.  The  Karhunen-Loeve  ex¬ 
pansion  along  with  Galerkin's  method  is  used  to  generate  a  reduced  order  model.  A  control  algorithm  is  achieved  by  applying 
Imear  quadratic  regulator  theory'  to  the  reduced  order  model. 

1  he  Karhunen-Loey  e  basis  functions  or  mode  shapes  must  first  be  found  to  identify  the  reduced  order  model.  Previous 
results  have  shov.  n  that  in  the  limit  as  the  structural  damping  approaches  zero  the  Karhunen-Eoeve  mode  shapes  and  natural 
mode  shapes  converge.  Numcncal  techniques  are  applied  to  evaluate  the  structural  damping  required  for  convergence.  Once 
the  Kariiuncn-Locve  mode  shapes  are  determmed.  the  reduced  order  control  model  is  applied  to  the  full  order  sy  stem.  I  he 
performance  of  various  Karhunen-Eoeve  models  is  compared  by  measuring  the  modal  energies  in  the  controlled  and  uncon¬ 
trolled  modes. 


Distnoulion  .A\ .T.is'y.kA  of  .'y'^sirav'. 


lla  Na:”0  'C  Rvspons 

-Jeff  B,  Burl 


I  :  .a?  r??or; 


oj  lodi'..d,.a: 


21  .Abslraci  Security  Classification 
□  DTic  user?  Unclassified 


22b  telephone  (include  .4rea  cede'  22c  Office  Svnibol 

(408)  646-2391)  62B1 


DU  t  f'/R  vi  1 4  s..’'-  MAR 


8.-  APR  edition  .Tiay  be  used  until  exhausted 
All  oiner  editions  arc  uosolete 


sccuriA  class. licuron 


I  nclassilicd 


1 


Approved  lor  public  release;  distribution  is  unlimited. 


Active  Damping  of  Vibration  in  Large  Space 
Structures  Using  a  Karhunen— Loevc  Reduced  Order  Model 

by 

Terence  M.  Grogan 
Lieutenant,  (,'nited  States  \a\T 
B.S..  Michigan  State  University,  19S0 

Submitted  in  partial  fulfillment  of  the 
requirements  for  the  degree  of 

MASTER  OF  SCIENCE  IN  ELECTRICAL  ENGINEERING 

from  the 

NAVAL  POSTGIUADUATE  SCHOOL 
March  1989 


Author; 


•Approved  by: 


John  P.  Powers.  Chairman. 

Department  of  Electrical  and  Computer  Engineering 


ABSTRACT 


Large  space  structures  are  difTicult  to  control  because  of  the  high  order  of  their 
mathematical  models.  The  high  order  mathematical  model  makes  the  use  of  a  reduced 
order  model  to  control  the  structure  desirable.  The  Karhunen-Loeve  expansion  along 
with  Galerkin  s  method  is  used  to  generate  a  reduced  order  model.  A  control  algorithm 
is  achieved  by  applying  linear  quadratic  regulator  theory  to  the  reduced  order  model. 

The  Karhunen-Loeve  basis  functions  or  mode  shapes  must  first  be  found  to  identify 
the  reduced  order  model.  Previous  results  have  shown  that  i^i  the  lintit  as  the  structural 
damping  approaches  zero  the  Karhunen-Loeve  mode  shapes  and  natural  mode  shapes 
converge.  Numerical  techniques  are  applied  to  evaluate  the  structural  damping  required 
for  convergence.  Once  the  Karhunen-Loeve  mode  shapes  are  determined,  the  reduced 
order  control  model  is  applied  to  the  full  order  system.  The  performance  of  various 
Karhunen-Loeve  models  is  compared  by  measuring  the  modal  energies  in  the  controlled 
and  uncontrolled  modes. 
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THESIS  DISCLAIMER 


The  reader  is  cautioned  that  computer  programs  developed  in  this  research  may  not  ♦ 

have  been  exercised  for  all  cases  of  interest.  While  eveiy  effort  has  been  made,  within 
the  time  available,  to  ensure  that  the  programs  are  free  of  computational  and  logic  er¬ 
rors.  they  cannot  be  considered  validated.  Any  application  of  these  programs  without 
additional  verification  is  at  the  risk  of  the  user. 
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I.  INTRODUCTION 


A.  LARGE  SPACE  STRUCTl  RES 

The  lightweight,  flexible  materials  used  to  construct  large  space  structures  (LSS), 
like  a  space  station,  are  lightly  damped  and  when  disturbed  will  vibrate  for  a  consider¬ 
able  amount  of  time.  This  prolonged  vibration  could  jeopardize  the  structural  integrity 
of  the  structure  or  disturb  experiments  on  the  LSS.  The  purpose  of  this  thesis  is  to  study 
the  effects  of  controlling  this  structural  vibration  with  a  reduced  order  model,  specifically 
a  Karhunen-Loeve  reduced  order  model. 

B.  PROBLEM  APPROACH 

1  he  solution  of  the  vibration  control  problem  requires  the  use  of  a  mathematical 
riiodel  that  describes  the  behavior  of  the  system  in  time.  The  space  structure  is  modeled 
as  a  combination  of  small  plates  of  unit  mass  connected  to  form  the  complete  structure. 
The  vibrational  motion  of  these  plates  can  be  modeled  as  a  set  of  coupled  damped  har¬ 
monic  oscillators.  Using  modal  analysis  the  mathematical  model  of  the  structure  can 
be  decoupled  to  yield  a  set  of  uncoupled  simultaneous  second  order  dilTerential 
equations. 

For  LSS  this  model  is  of  very  high  order,  making  control  design  difficult.  It  is 
therefore  desirous  to  use  a  reduced  order  model  (ROM)  to  control  the  structure.  A 
control  system  is  designed  for  the  LSS  using  the  Karhunen-Loeve  (reduced  order)  model. 
The  control  system  is  then  applied  to  the  LSS  after  it  has  been  disturbed  by  an  impulse. 
The  performance  of  this  system  and  a  control  system  based  on  a  modal  model,  truncated 
by  frequency,  are  compared. 

The  LSS  used  in  this  thesis  is  a  dual  keel  space  station,  see  Figure  1  on  page  2.  A 
mathematical  model  of  this  space  station  is  provided  courtesy  of  McDonnel  Douglas 
.Astronautics.  .A  computer  simulation  of  this  space  station  is  used  to  determine  the  ef¬ 
fectiveness  of  the  control  system.  Karhunen-Loeve  (reduced  order)  models  of  increasing 
size  are  simulated  and  the  modal  energies  are  calculated  and  used  to  determine  the  rela¬ 
tive  cflectiveness  of  the  control  models.  These  results  will  then  he  compared  to  results 
obtained  for  the  modal  model  truncated  by  frequency. 


Figure  1.  Representation  of  a  Dual  Keel  Space  Station  from  [Ref.  1:  p.  1] 


C.  ORG.4iNIZATION 

•  Chapter  II  is  a  statement  of  the  mathematical  model  of  the  space  station. 

•  Chapter  III  discusses  the  Karhunen— Loeve  expansion,  the  relationship  between 
the  Karhunen-Loeve  mode  shapes  and  the  natural  mode  shapes  and  how  a  re¬ 
duced  order  model  can  be  synthesized  using  the  Karhunen-Loeve  expansion  and 
Galerkin's  method. 

•  Chapter  IV  discusses  the  determination  of  the  Karhunen-Loeve  modes. 

•  Chapter  V  describes  the  control  solution  and  presents  the  control  simulation  re¬ 
sults. 

•  Chapter  VI  discusses  the  conclusions  based  on  the  simulation  results  and  a  re¬ 
commended  design  procedure. 


II.  THE  MATHEMATICAL  MODEL 

A.  INTRODCCTION 

Large  space  structures  are  flexible,  lightly  damped  structures  that  can  vibrate  for  a 
considerable  amount  of  time  when  disturbed  by  an  external  force.  To  effectively  control 
such  a  system,  a  mathematical  model  describing  the  evolution  of  that  system  in  time  is 
required. 

B.  THE  MODAL  MODEL 

The  space  station  is  modeled  as  a  combination  of  small  plates  of  unit  mass  con¬ 
nected  to  form  the  complete  structure.  This  model  can  be  visualized  as  a  system  of 
masses  connected  by  springs  and  dashpots  (the  springs  representing  the  stiffness  factor 
and  the  dashpots  representing  the  damping  factor).  The  displacement  of  the  masses  can 
be  described  by  a  second  order  matrix  differential  equation  of  motion: 

K.\ii)  +  KxiO  =  F(i)  { 1 ) 

where 

•  \  is  the  generalized  coordinate  vector 

•  M  is  the  diagonal  system  mass  matrix 

•  K  is  the  structural  damping  term 

•  d  is  the  damping  constant 

•  (:>.■  is  the  frequency  of  oscillation  of  the  system 

•  K  is  the  symmetric  system  stiffness  matrix 

•  f(t)  is  the  system  forcing  function 

This  equation  represents  a  system  of  simultaneous,  second  order  differential  equations 
that  are  coupled  by  the  K  matrix. [Ref  1:  p.  3] 

1  his  equation  can  be  uncoupled  and  the  system  represented  by  a  set  of  independent 
second  order  differential  equations.  This  is  done  through  the  process  of  modal  analysis 
which  is  outlined  ’  [Refs.  2.  3].  The  resulting  modal  model  consists  of  a  set  of  inde¬ 
pendent  second  order  differential  equations: 

[  i/  -f  Q  4-  j/  =  T]  ( 2 ) 
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\\here 


•  t]  is  the  coordinate  vector  or  modal  amplitude  vector 

•  Cl-  =  diag^fc)^,,  ...  ,  aj,J 

•  A’''=  C-v,  .Vi  ...  .V,]  the  transpose  of  the  modal  matrix  or  mode  shape  vector 

•  F  is  the  torquing  force  applied  at  a  point 

Next  a  discrete-time  equation  describing  the  motion  of  the  space  station  in  terms  of  its 
natural  modes  of  vibration  is  developed.  The  discrete-time  state  equation  for  the  ith 
equation  of  motion  is; 

Z;(/cr+  i)  =  <t),(7iz,(^r)-hr,(7ixf[n^7i  +  n'(A-7^]  (3) 


where 

•  Z,  is  a  vector  of  the  ith  modal  amplitude  and 
the  ith  modal  velocity 

•  O  is  the  ith  state  transition  matrix 

•  r,  is  the  ith  input  vector 

•  Xf  is  the  transpose  of  the  ith  mode  shape  vector 

•  F  is  the  control  torque  force  vector  applied  at  a  point 

•  T  is  the  sampling  time 

•  k  is  the  time  inde.x 

•  W  is  the  disturbance  input 

This  equation  is  used  for  computer  simulation  of  the  space  station  and  control 
solution. [Ref  1:  p.  4] 
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III.  APPLICATION  OF  THE  KARHUNEN-LOEVE  EXPANSION  TO  THE 
REDUCED  ORDER  CONTROL  OF  LARGE  SPACE  STRUCTURES 


A.  INTRODUCTION 

Large  flexible  structures,  such  as  a  space  station,  as  a  class  of  distributed  parameter 
systems  (DPS)  require  a  finite  dimensional  model  for  control  design.  This  model  may 
be  achieved  by  approximating  the  state  of  the  LSS  using  the  Karhunen-Loeve  (KL)  ex¬ 
pansion.  The  expansion  is  truncated  to  provide  the  finite  dimensional  approximation 
of  the  state  for  control  design.  The  KL  model  that  results  describes  the  evolution  of  the 
approximated  state  of  the  structure. 

The  natural  mode  shapes  are  normally  used  for  modeling  and  control  of  flexible 
structures.  The  relationship  between  the  natural  mode  shapes  and  the  KL  mode  shapes 
is  described  in  this  chapter  as  well  as  the  use  of  Galerkin's  method.  Galerkin's  method 
is  used  to  generate  a  reduced  order  model  (ROM),  using  both  the  natural  mode  shapes 
and  the  KL  expansion. 

B.  THE  KARHUNEN-LOEVE  EXPANSION 

The  purpose  of  the  KL  expansion  is  stated  by  Stark  and  Woods  [Ref  4;  p.  322]. 
“The  idea  [of  the  KL  expansion]  is  to  decompose  a  general  second-order  random  proc¬ 
ess  into  an  orthonormal  expansion  whose  coelTicients  are  uncorrelated  random 
variables.''  The  state  of  a  large  space  structure  (LSS),  y(x),  can  be  modeled  as  random 
process  since  it  depends  on  random  excitations,  i.e.,  noise  from  onboard  machinery  and 
actuators.  The  second  order  moments  of  the  LSS  are  proportional  to  the  physical  en¬ 
ergy  and  therefore  exist.  The  LSS  can  therefore  be  approximated  using  the  KL  expan¬ 
sion  [Ref  5:  p.  12].  This  is  done  by  projecting  the  random  process  onto  an  orthonormal 
basis  and  truncating  to  N  terms.  The  value  chosen  for  N  is  a  matter  of  “engineering 
judgment".  [Ref  5;  p.  13] 

The  selection  of  the  orthonormal  basis  is  made  by  solving  the  KL  eigenequation: 

<  Ryy{.x,z),  (t>iiz)  >  =  (4) 


where 

•  0,(.r)  is  referred  to  as  an  eigenfunction  (or  the  KL  mode  shapes) 

•  /,  is  the  eigenvalue  and  is  a  measure  of  the  excitation  of  the  ith  basis  function 
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•  /?^,(.r.z)  is  the  correlation  function  of  y(x):  Ryy{x,z)  = 

•  <  • ,  •  >  is  an  inner  product;  <  a(z),  b(z)  >  =  M(z)  b(z)  dz 

•  M(z)  is  the  mass  density  of  the  structure 

•  Q  is  the  spatial  extent  of  the  structure 

The  state  of  the  LSS  is  approximated  b\': 


(5) 


where 

•  is  an  approximation  to  the  state  of  the  space  structure 

•  C,  is  a  set  of  coordinates  found  by  c,  =  <  4>,{x),y{x)  > 

•  4>,{x)  is  the  ith  basis  function  or  KL  mode  shapes 

The  expansion  is  truncated  by  keeping  the  eigenfunctions  (KL  mode  shapes)  associated 
with  the  N  largest  eigenvalues.  [Ref  5:  p.  13] 

The  KL  expansion  yields  the  best  approximation  to  the  random  process,  i.e.,  mini¬ 
mizes  the  expected  value  of  the  norm  of  the  error,  of  any  orthogonal  e,xpansion.  The 
appro.xin'ation  error  is  defined  as: 

OO 

c(.r)=y(.Y)->aW=  ^  W  (6) 

/=.v+i 

For  a  proof  of  the  optimality  of  the  KL  expansion  see  [Ref  6:  p.  11]  or  [Ref  5:  p.  15]. 

C.  RELATIONSHIP  BETWEEN  THE  KARHUNEN-LOEVE  MODE  SHAPES  AND 
THE  NATURAL  MODE  SHAPES 

The  relationship  between  the  KL  mode  shapes  or  KL  basis  functions  and  the  natural 
mode  shapes  is  taken  from  Burl  [Ref  5:  p.  13].  The  KL  mode  shapes  of  a  structure  are 
related  to  the  natural  mode  shapes  of  that  structure  by  a  linear  transformation  repres¬ 
enting  a  change  of  basis  which  can  be  written: 


(7) 


6 


where 


0 

V2M 

0 

(8) 

0 

r/M) 

0 

>/2(-^) 

0 

- 

T  r  1 
q  =|c,. 

2  3 

q  q 

1 

(9) 

and  {ri,{x)}  is  the  set  of  natural  mode  shapes.  The  state  of  a  structure  consists  of  a 
generalized  position  and  velocity  which  can  be  expanded  in  terms  of  the  natural  mode 
shapes: 


r(.r,/)  =  y  a//) 


7=1 


vp) 

0 


+  ay(t) 


0 

Vp) 


(10) 


where 


a^(r)  =  [a,(/)  a,(/)  a^it)  otjC/)  a3(/)  ...  ]  (11) 

are  the  coordinates  and  velocities  of  the  natural  mode  shapes.  The  vectors,  c,.  can  be 
found  by  solving  the  equation: 

£[«(/)  a^(/)]Qc,  =  ;.,c,  (12) 


where 


Q 


=  diag 


COy  0 

r  7  rti 

cog  0 

0 

1 

0 

_ 1 

- 1 

0 

_ 1 

- 1 

0 

_ 1 

(13) 


Equations  10  and  11  give  the  KL  mode  shapes  in  terms  of  the  natural  mode  shapes. 
Equation  12  is  an  infinite  dimensional  eigenvalue  problem  that  can  be  solved  practically 
by  truncating  it  to  the  most  significant  terms. 

D.  REDUCED  ORDER  MODELING  USING  THE  KARHUNEN-LOEVE  MODES 

Galerkin's  method  is  used  to  produce  reduced  order  state  equations.  This  can  be 
done  with  either  the  KL  modes  or  the  natural  modes. 
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The  discrete-time  state  equation  of  a  distributed  parameter  system  is: 


y{.x,k)  =  -  1)  -h  GAk) 


(!4) 


where 

•  y(x,k)  is  the  state 

•  flk)  is  the  input 

•  k  is  the  time  index 

•  F  is  an  operator  on  the  state  space 

•  G  is  an  operator  from  to  the  state  space 

•  For  each  k,  y{x,k)  e  the  state  space,  f{k)  e 

A  finite  dimensional  approximation  to  this  equation  can  be  obtained  using  Galerkin's 
method 

P„y{x,k)  =  {P„FP„)yix,k  -  1)  -h  G  u{k)  (15) 

where  P„  can  be  written  in  terms  of  a  basis  1  '  =  1,2, ... ,  «!] 

n 

Pni‘)  =  'yjViix)<  Vi{x),  •>  (16) 

i=l 

The  KL  mode  shapes  (basis  functions)  can  be  used  in  the  above  equations  and  the  KL 
model  results.  The  natural  mode  shapes  can  be  used  in  the  above  equations,  which  is 
equivalent  to  truncating  the  modal  equations,  producing  the  modal  model. [Ref  5  :  p. 
131 

E.  SUMMARY 

The  KL  expansion  can  be  used  to  approximate  the  state  of  a  LSS.  This  is  done  by 
projection  of  the  state,  which  is  modeled  as  a  random  process,  onto  an  orthonormal 
basis  function.  The  basis  function  can  be  found  by  solving  the  KL  eigenequation.  These 
KL  basis  functions  will  yield  the  best  approximation  to  the  state  of  the  LSS.  Then  using 
Galerkin's  method  a  reduced  order  model  of  the  LSS  is  produced.  This  reduced  order 
model  is  used  to  generate  a  control  which  is  applied  to  the  entire  system. 
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IV.  DETERMINATION  OF  THE  KARHUNEN-LOEVE  BASIS 

FUNCTIONS 


A.  INTRODUCTION 

The  KL  basis  functions  (or  mode  shapes)  can  be  determined  by  solving  the 
eigenequation  12.  This  is  a  tedious  and  laborious  processes.  It  can  be  shown  that  for 
lightly  damped  structures  the  KL  mode  shapes  can  be  determined  from  the  open  loop 
response  by  ordering  the  natural  mode  shapes  in  order  of  decreasing  modal  energies. 

The  KL  mode  shapes  were  calculated  using  the  KL  mode  program  in  Appendix  B. 
These  calculations  were  compared  to  the  mode  shapes  selected  by  observing  the  open 
loop  response  to  verify  that  this  is  a  valid  method  of  determining  the  KL  mode  shapes. 

The  value  of  the  structure's  damping  factor  used  in  the  program  was  increased  until 
the  KL  mode  shapes  no  longer  were  the  same  as  those  determined  from  the  open  loop 
response.  This  was  done  to  determine  how  large  the  structure's  damping  factor  could 
be  and  still  have  the  KL  mode  shapes  converge  to  the  natural  mode  shapes. 

B.  NUMERICAL  DETERMINATION  OF  THE  KARHUNEN-LOEVE  MODE 
SHAPES 

The  KL  mode  shapes  are  found  by  determining  the  eigenvalues  and  eigenvectors 
(equation  12)  of  the  covariance  matrix  £[«,(/)  This  is  an  infinite  dimensional 

matrix  which  is  truncated  to  a  finite  dimensional  square  matrix.  There  are  three  terms 
in  E\iy.ii)  where  a[i)  is  defined  by  equation  11.  The  first  of  these  terms  is  com¬ 

puted.  for  a  white  noise  input 

E[«,(t)«/t):=  [^(Tih/Tior  (17) 

''0 

where  h,(t)  is  the  impulse  response  for  a,(t)  [Ref  6:  p.  66].  The  impulse  response  is  given 
by: 


h;(f)=' 


''2/ 


CO, 


oi\ 


1 


1  -  y  t) 


(18) 


where 
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•  bj,  is  the  modal  slope  in  the  x  direction 

•  y  is  the  damping  coefficient;  y  =  -^  is  a  constant 

•  is  the  natural  frequency 

Performing  the  integration  in  equation  17  produces  the  result: 


K,K,[  2[a^  +  c^]] 


where 


•  a  =  y  w,,  +  y  w,, 

•  b  =  OJ,,  ^  1  -  y*  -  1  -  y* 

•  c  =  1  -  y^  +  1  -  y* 

The  other  terms  in  E[o!(t)  aTt)3  can  be  found  in  a  similar  way.  They  are; 


a,(t)]  =  K,  K,  L,. 


.  2[a'  +  b=]  2[a=  +  c=] 


ECo:,(I)  a,(l)]  -  K,  K 


2[a^  +  b^]  2[a^  +  c^] 


L,L,a  —  LjM.b  —  L,M,b  +  MjM;a 


2[a'  +  b'] 

EjLya  -  L,MyC  —  M,Mya 

2[777] 


where 

•  L,,^  =  T  -  y^ 

•  M,,,  =  y  w,,., 

Equations  19,  20,  and  21  specify  the  eigenvalue  problem  equation  12. 

The  eigensolution  specifies  the  transformation  from  the  natural  mode  shapes  to  the 
KL  mode  shapes.  This  finite  dimensional  eigenvalue  problem  is  solved  and  the  KL  mode 
shapes  are  approximated  as  a  linear  combination  of  the  first  fifty  (flexible)  natural  mode 


10 


shapes.  The  program  in  Appendix  B  computes  the  KL  mode  shapes  by  solving  the 
eigenvalue  problem  equation  12. 

C.  EMPIRICAL  DETERMINATION  OF  THE  KARHUNEN-LOEVE  MODE 
SHAPES 

The  numerical  computation  of  the  KL  mode  shapes  is  difficult.  Solving  the 
eigenvalue  problem  requires  the  deternrination  of  the  impulse  response  and  solving  for 
all  the  terms  in  the  covariance  matri.x  Ena,(t)  a,,(t)].  It  turns  out  that,  in  the  limit  as  the 
damping  factor  of  the  structure  approaches  zero,  the  KL  mode  shapes  converge  to  the 
natural  mode  shapes  [Ref  7]. 

If  the  modal  energies  for  the  open  loop  response  are  computed  and  the  mode  shapes 
ordered  by  decreasing  modal  energy,  these  mode  shapes  should  be  the  same  as  those 
computed  numerically.  The  output  of  the  system,  y,  is  defined: 


OJj 

0 

0 

0 

”  y'7 ' 

0 

1 

0 

0 

'^7  ' 

• 

• 

• 

•  • 

• 

• 

• 

• 

♦ 

•  • 

• 

• 

« 

• 

• 

•  • 

• 

• 

356 

0 

0 

0 

•>^56 

0 

0 

0 

1 

where  y,  e  d?-  and 

£[ith  .Modal  Energy]  =  £[  Wy^  |p]. 
The  energy,  given  in  equation  23  can  be  written  [Ref  S  ]: 


(22) 


(23) 


roo 

£IILvvll']=  I  WhyMfoi  (24) 

•'0 

where  h^,  is  the  response  ofy,  due  to  an  impulse  applied  at  the  disturbance  input  {node 
69  or  55).  Equation  24  is  evaluated  using  computer  simulation.  The  mode  shapes  found 
are  compared  to  those  determined  numerically  for  noise  input  at  node  69  or  node  55,  see 
Table  1  on  page  12. 
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Table  1  shows  that  the  mode  shapes  are  the  same.  Therefore,  the  KL  mode  shapes  can 
be  determined  from  the  open  loop  response  when  the  structural  damping  factor  is  suffi- 
ciently  small.  In  the  case  of  the  space  station  simulated  in  this  thesis  a  value  of  0.001 
was  used.  The  next  question  is.  how  big  can  the  damping  factor  of  the  structure  be  be¬ 
fore  the  KL  mode  shapes  fail  to  converge  to  the  natural  mode  shapes? 

1.  Required  Magnitude  of  the  Damping  Factor  for  Convergence 

When  the  structural  damping  factor  is  small  enough,  the  eigenvectors  solved  for 
in  equation  12  approximate  the  natural  basis 

e|  =  Cl00...  ]  (25) 


or 


eT  =  [0  1  00  ...  ] 


(26) 


etc.;  one  element  is  unity  and  the  other  elements  are  zero.  When  these  vectors  are  sub¬ 
stituted  into  equation  7  it  is  easy  to  see  that  the  KL  mode  shapes  are  the  natural  mode 
shapes.  If  the  damping  factor  is  increased,  at  some  value  the  eigenvectors  have  more 
than  one  non-zero  element  and  the  KL  mode  shapes  become  a  linear  combination  of 
natural  mode  shapes. 

The  damping  factor  was  increased  successively  by  a  factor  of  two  from  a  starting 
value  of  0.0005.  The  norm  of  the  error  (or  error  norm)  is  used  as  a  means  of  measuring 
how  closely  the  eigenvectors  obtained  from  the  KL  mode  program  approach  the  ideal 
of  equation  25.  The  error  norm  is  defined  as  follows: 


II  £  lljj  =  minimumll  —  ej  ||  = 


(27) 


where 

4  =  [Xl  ^2  X3  -  Xj  ...  X50] 


(28) 


and  for  simplicity  x,  is  assumed  to  be  positive  and  the  largest  component  ofc,.  Equation 
27  can  be  simplified  to 
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Since  the  eigenvectors  are  normalized. 


50 

=  1  (30) 

/=i 

and  substituting  equation  30  into  equation  29  yields 

|U||,  =  ^2(l-xp  (31) 

a.  Results 

The  norm  of  the  error,  t,  and  the  natural  mode  shape  sequence  obtained 
from  the  KL  mode  program  were  used  to  determine  how  large  the  damping  factor  could 
be  and  still  have  the  KL  mode  shapes  and  natural  mode  shapes  converge.  The  error 
norm  is  presented  graphically  for  each  KL  mode;  all  cases  are  for  disturbance  torques 
due  to  actuator  noise.  The  cases  presented  are: 

•  Tigure  2  on  page  16:  for  a  damping  factor  of  0.0005  all  but  two  KL  modes  have 
an  error  norm  below  0.2.  The  natural  mode  sequence  is  as  in  Table  1  on  page  12. 

•  Figure  3  on  page  17;  for  a  damping  factor  of  0.001  nearly  74  percent  of  the  KL 
modes  have  an  error  norm  of  0.2  or  below,  the  natural  mode  sequence  is  presented 
in  Table  1  on  page  12.  The  natural  mode  shapes  are  a  good  approximation  of  the 
KL  mode  shapes. 

•  Figure  4  on  page  18;  for  a  damping  factor  of  0.002  only  46  percent  of  the  KL 
modes  have  an  error  norm  at  or  below  0.2.  and  the  natural  mode  sequence  no 
longer  conforms  to  that  in  Table  1  on  page  12,  i.e.,  the  KL  modes  and  natural 
modes  are  not  converging. 

•  Figure  5  on  page  19;  for  a  damping  factor  of  0.005  only  28  percent  of  the  KL 
modes  have  an  error  norm  of  0.2  or  less  and  the  KL  modes  and  the  natural  m.odes 
are  more  divergent. 

b.  Conclusion 

If  the  damping  factor  is  greater  than  0.001  the  natural  mode  shapes  are  not 
a  good  approximation  of  the  KL  mode  shapes.  The  convergence  criteria  is  that  the 
norm  of  the  error,  for  50  percent  of  the  KL  modes  or  greater,  is  more  than  0.2.  The 
value  for  the  damping  factor  used  in  the  space  station  simulations  as  noted  in  previous 
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chapters  was  0.001.  The  a.ssumption  that  the  KL  mode  shapes  and  natural  mode  shapes 
converge  is  valid  for  this  damping  factor  and  determining  the  KL  modes  from  the  modal 
energies  is  a  good  approximation. 
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ERROR  NORM  VS  DAMPING  FACTOR 


o 


MODE 


Figure  2.  Norm  of  the  error  values  for  the  KL  modes,  excitation  at  node  69  , 
damping  factor  =  0.0005 
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ERROR  NORM  VS  DAMRING  FACTOR 
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ERROR  NORM  VS  DAMPING  FACTOR 
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Figure  4.  Norm  of  tlie  error  values  for  the  KL  modes,  excitation  at  node  69  , 
damping  factor  =  0.002 
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Norm  of  the  error  values  for  the  KL  modes,  excitation  at  node  69 
damping  factor  =  0.005 


V.  CONTROL  SOLUTION 


A.  INTRODUCTION 

The  response  of  the  space  station  to  disturbance  and  control  inputs  is  simulated 
using  the  modal  model  discussed  in  Chapter  II.  The  Fortran  program  that  simulates  this 
model  was  written  by  Preston  (Ref.  I;  p.  47],  and  is  used  in  this  thesis  with  minor  mod¬ 
ification,  (see  .Appendi.x  A).  The  objective  of  the  simulation  is  to  determine  the  system 
response  to  disturbances  applied  at  modes  69  and  55  (see  Figure  6)  on  the  structure, 
using  increasing  sizes  of  KL  models  (i.e.,  5,  10,  20  modes).  The  response  of  the  system 
is  depicted  graphically  by  displaying  the  energy  in  each  mode  in  english  units  of  inch- 
pounds  (in— lbs). 


B.  FULL  ORDER  SPACE  STATION  MODEL 


The  plant  being  controlled  in  this  thesis,  as  stated  in  Chapter  1,  is  a  preliminary’ 
version  of  the  NASA  dual  keel  space  station.  NASTR.AN  was  used  to  generate  the  first 
fifty  ficxible  modes  (starting  with  mode  seven).  The  modal  model  that  forms  the  full 
order  model  is  composed  of  these  flexible  modes; 
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modal  amplitude 
modal  velocity 


(33) 


and 

•  0),  is  the  natural  frequency  of  the  ith  mode 

•  d  =  0.001  is  the  damping  coefficient 

•  u  6  .^Ms  the  control  input  torques  from  three  orthogonally  oriented  control  mo¬ 
ment  gyros 

•  ;v  6  ^Ms  a  random,  white  noise  disturbance  input  {three  orthogonal  disturbance 
torques) 

•  4>,{p)  e  is  the  modal  slopes  at  p 

•  is  the  location  of  the  control  moment  gyros  (node  69  on  Figure  6  on  page  20) 

•  />„  is  the  location  of  the  disturbance  input  (either  the  control  moment  gyro  location, 
node  69.  or  the  alpha  joint,  node  55,  as  seen  in  Figure  6  on  page  20). 
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This  model  is  used  as  the  full  order  model  in  the  simulation  of  the  space  station.  How¬ 
ever.  it  is  itself  a  reduced  order  model  since  the  “infinite"  number  of  vibrational  modes 
of  the  structure  are  truncated  to  fifty. 

It  is  assumed  that  perfect  modal  amplitude  and  modal  velocity  information  is  avail¬ 
able.  This  assumption  of  perfect  sensor  information  is  not  realistic,  but  isolates  the  ef¬ 
fect  of  modal  truncation  on  control  algorithm  synthesis. 

C.  THE  REDUCED  ORDER  MODEL 

The  reduced  order  model  of  the  space  station  is  obtained  by  truncating  modes.  The 
criteria  for  truncating  modes  is  either  the  modal  frequency  or  the  Karhunen-Loeve  or¬ 
dering  (the  KL  ordering  found  in  the  previous  chapter  to  be  equivalent  to  an  ordering 
based  on  the  energy  measured  during  open  loop  excitation). 

D.  PERFORMANCE  FUNCTION  AND  OPTIMAL  CONTROL 

The  full-order  model  and  a  reduced  order  model  have  been  established.  Next,  a 
mathematical  expression  of  system  performance  is  required.  The  performance  function. 
J.  stated  here  is  developed  in  detail  in  [Ref  1;  p.  10  ].  J  consists  of  the  total  energy  in 
the  modeled  modes  plus  a  control  energy  term: 

0  g' 

.  ^  0  r  -^7  ■ 

’  •  •  r 

*  •  •  +  /'w  w  (.>4) 

•  2  Y 

0  1 

where 

•  imax  is  the  number  of  modes  in  the  reduced  order  model 

•  r  =  10->7  is  a  weighting  coefficient  on  the  control  energy  term  selected  to  yield  time 
constants  on  controlled  modes  of  appro.ximately  30  seconds. 

Vibration  damping  is  achieved  by  application  of  steady  state,  linear  quadratic  regulator 
theory  [Ref  9]  to  the  KL  (reduced  order)  model.  The  control  torque  vector,  u(k),  is  the 
product  of  an  optimal  gain  matrix.  L,  and  the  time  varying  state  matrix,  Z  (equation  3). 
The  L  matrix  is  found  by  solution  of  the  Ricatti  equations  to  minimize  the  performance 
function.  J,  (Ref  1:  p.  12]. 


E.  SIMULATION  RESULTS 


The  space  station  was  simulated  and  the  modal  energies  calculated  for  disturbance 
inputs  applied  at  nodes  69  and  55  as  shown  in  Figure  6  on  page  20.  The  contrc^  system 
consists  of  KL  (reduced  order)  models  of  five,  ten  or  twenty  modes.  The  data  is  pre¬ 
sented  graphically  as  energy  in  each  of  the  modes.  The  open  loop  response  is  included 
for  comparison  except  where  the  difference  in  scale  precludes  it.  The  closed  loop  re¬ 
sponse  for  a  reduced  order  model  generated  by  modal  truncation  is  included  for  com¬ 
parison  with  the  KL  model.  The  cases  presented  are; 

•  Open  loop  (no  control)  response  which  identifies  the  reduced  order  model  to  be 
used.  See  Figure  7  and  Figure  S. 

•  Closed  loop  response  for  a  reduced  order  model  generated  by  modal  truncation 
with  the  modes  ordered  by  natural  frequency.  See  Figure  9  to  Figure  14. 

•  Closed  loop  response  for  a  Karhunen-Loeve  (reduced  order)  model.  The  KL  modes 
were  determined  from  the  open  loop  responses.  These  results  are  shown  in 
Figure  15  through  Figure  19. 

•  Closed  loop  response,  with  a  KL  model  based  on  the  open  loop  response,  for  the 
system  excited  by  the  control  actuators.  These  results  are  shown  in  Figure  20 
through  Figure  22. 

These  results  are  discussed  in  detail  in  Chapter  VI. 
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Figure  7.  Open  Loop  Response,  excitation  at  node  55. 
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Figure  8.  Open  Loop  Response, excitation  at  node  69 
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9.  Five  Controlled  modes  selected  by  frequency,  excitation  at  node  5 
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Figure  10.  Ten  controlled  modes  selected  by  frequency,  excitation  at  node  55. 
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Figure  11.  Twenty  controlled  modes  selected  by  frequency, e.xcitation  at  node  55. 
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Figure  12.  Five  controlled  modes  selected  by  frequency,  excitation  at  node  69. 
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Figure  13.  Ten  controlled  modes  selected  by  frequency,  e.xcitation  at  node  69. 
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Figure  14.  Tnenty  controlled  modes  selected  by  frequency,  e.xcitation  at  node  69. 
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Figure  15.  Five  controlled  modes  selected  by  modal  energy,  e.xcitation  at  node  55. 
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Figure  16.  Ten  controlled  modes  selected  by  modal  energy,  e.xcitation  at  node  55. 
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Figure  17.  Twenty  controlled  modes  selected  by  modal  energy,  excitation  at  node 
55. 
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Figure  18.  Five  and  ten  controlled  modes  selected  by  modal  energy,  excitation  at 


node  69, 
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ENERGY  (IN-LBS)  *10"' 


Figure  20.  Five  controlled  modes  selected  by  modal  energy  due  to  excitation  at  node 
69,  actual  excitation  at  node  55. 
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Figure  21.  Ten  controlled  modes  selected  by  modal  energy  due  to  excitation  at  node 
69,  actual  excitation  at  node  55. 
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Figure  22.  Twenty  controlled  modes  selected  by  modal  energy  due  to  excitation  at 
node  69,  actual  excitation  at  node  55. 
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VI.  CONCLUSIONS  AND  RECOMMENDATIONS 


A.  CONCLUSIONS 

Results  supported  theoretical  claim  that  the  KL  mode  shapes  converge  to  the  na¬ 
tural  mode  shapes.  For  the  space  station  example,  a  damping  factor  of  0.001  was  re¬ 
quired  to  make  the  approximation  based  on  convergence  good. 

Preston  discovered,  see  [Ref  1:  p.  46],  that  truncating  the  modes  based  on  natural 
frequency  yielded  poor  results.  Modes  with  large  coupling  to  the  control  system  are  not 
modeled  when  natural  frequency  is  the  method  of  truncation.  Instead  modes  with  verv' 
little  coupling  to  the  control  system  are  modeled  and  in  attempting  to  control  these 
modes  ver\'  large  control  torques  are  generated  causing  large  excitations  in  the  strongly 
coupled  but  unmodeled  modes  as  seen  in  Figures  9  through  14. 

The  excitation  of  each  mode  depends  on  two  factors:  the  natural  frequency  which 
determines  the  damping  and  the  amplitude  of  the  mode  shape  which  determl’^es  how 
much  excitation  is  received  by  each  mode.  The  amplitude  of  the  mode  shape  is  the 
dominate  factor  which  is  evident  from  the  open  loop  responses  see  Figures  7  and  8. 

Truncating  based  on  the  open  loop  response,  i.e.,  the  Karhunen-Loeve  model  also 
yields  poor  results  except  when  it  is  developed  from  the  open  loop  response  for  the  case 
of  the  disturbance  torques  being  due  to  actuator  noise,  see  Figure  8  on  page  25.  When 
truncating  modes  based  on  the  open  loop  response,  the  coupling  from  the  noise  input 
is  the  dominate  factor  in  selecting  the  modes  to  include  in  the  model.  In  this  case  modes 
with  large  open  loop  excitation  are  controlled,  but  the  control  coupling  still  dontinates 
yielding  large  excitation  in  the  unmodeled  modes  (see  Figures  15  through  17).  When  the 
disturbance  torques  are  due  to  actuator  noise,  the  modes  with  the  largest  open  loop 
excitation  are  also  the  modes  with  the  largest  control  coupling  and  the  control  system 
works  well  (see  Figures  18  and  19). 

The  KL  model  based  on  the  open  loop  response  for  disturbance  torques  due  to 
actuator  noise  works  well  when  applied  to  the  space  station  with  a  disturbance  torque 
applied  at  another  location.  As  seen  in  Figures  20  through  22,  the  control  system  drives 
the  modeled  modes  close  to  zero  without  exciting  other  modes  to  a  significant  degree. 
The  problem  with  this  configuration  is  that  modes  with  large  open  loop  excitation  are 
not  necessarily  modeled  and  may  be  left  excited  by  the  control  system. 


B.  DESIGN  PROCEDURE 

This  design  procedure  is  based  on  the  above  conclusions  and  is  recommended  for 
controlling  vibration  in  LSS.  First,  the  reduced  order  control  model  should  be  based  on 
the  KL  mode  shapes  determined  from  the  open  loop  response  for  disturbance  torques 
due  to  actuator  noise.  Next,  check  the  open  loop  response  for  disturbance  torques  at 
locations  on  the  structure  where  noise  inputs  are  most  likely,  e.g.,  the  alpha  joint. 
Modes  strongly  coupled  to  the  disturbance  input  must  be  included  in  the  KL  model. 
This  will  determine  the  size  of  the  KL  (reduced  order)  model.  If,  after  all  this  is  done, 
the  results  are  still  unacceptable,  then  adding  additional  control  elements  should  be 
considered.  The  placement  of  the  additional  control  element(s)  is  determined  by  the 
node(s)  with  the  largest  modal  amplitude(s)  for  the  modes  that  remain  excited  by  the 
initial  control  system.  The  process  is  repeated  until  an  acceptable  control  is  achieved. 
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APPENDIX  A.  SPACE  STRUCTURE  SIMULATION  PROGRAM 


This  program  simulates  the  dual  keel  space  station  described  in  Chapter  1  by  im¬ 
plementing  the  model  described  in  Chapter  IV.  The  control  described  in  Chapter  IV  is 
also  simulated.  For  a  detailed  explanation  of  this  program  and  its  development  see  [Ref. 
!]• 


*********************************************************** 


Space  Structure  Simulation  Program 


*****  By  William  J.  Preston  *v-*** 


y«.u*.u«.uy« 

**Vf*Vf****' 


VARIABLE  DECLARATIONS 


■*• 


EXTERNAL  EXCMS.RICDSD 
CHARACTER*6  NAM 

CHARACTER* 1  AGAIN , CORECT , RAGAIN 

INTEGER  NODE , MODE , KQ , EMODE , SMODE 

INTEGER  CT , CF , KA , LOOP , PRNT , MODAL , V , COUNT , PRNTG 

I.NTEGER  NF ,  NG ,  NH ,  NZ ,  I ,  K ,  M , CTAD J 

INTEGER  IPVS(100),ITYPE(200) 


REAL  T0TCST,RM0DEN(7:  100) 

REAL*8  COSWIT ,  S  IN^W IT ,  COST ,  CNTCST ,  ENERGY ,  RM 


REAL  RT0TAL,RM0DEN(7:  100) 

REAL*8  PHII(2,2,100),GAMMA(2,100),EGT,GMA,WN,W1,X1T,X2T 
REAL*8  PHI(188,188),B(188,3),BN(188,3),R(3,3),RR(3,3) 

REAL*8  RINVC3,3) ,RRINV(3,3),X1(7: 100),X2(7: 100) ,M0DEN( 7: 100) 
REAL*8  COSW IT , S INW IT , COST , CNTCST , ENERGY , TOTCST , RM 
REAL*8  TCX , TCY , TCZ , DAMP , SAMPT , PI , SUMl , SUM2 , SUM3 , SUMC 
REAL*8  TNX ,TNY ,TNZ , IMPX , IMPY , IMPZ , IMPLSX, IMPLSY , IMPLSZ 
REAL  LAMA( 100 ) , UGVEX( 100,3) ,RNODE , RHODE , MIN .TIME , SAMPTM 
REAL  UG69(100,3),UG23( 100, 3), UG55( 100,3) 

REAL*8  H( 100, 100) ,G( 100,100) ,L( 3, 100) ,BT( 3 , 100) 

REAL*8  Z(200,200),W(200,200),ER(200),F(100,100),EI(200) 
REAL*8  SCALE( 200) ,TEMP( 100,3) ,TEMP1( 3,100) ,WORK( 100) 
**■!>*************■>'?****************■************************** 


*********************************************************** 


*****  VARIABLE  DEFINITIONS  ***** 

*VfiV‘>VVfyf’>VTVyrTVff’}V^V'f')VVfiVVr*Vcyf**fp'}VyrVfTVTV^V^V>VVcyf'j^>V>VV'r*'A’^>V’>VTV’»VVp*VfyrycTVVp’>VTVVr**’>V‘>V 
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LAMA  =  VECTOR  OF  THE  SQUARE  OF  THE  NATURAL  FREQUENCIES 

UGVEX  =  MODE  POSITONS  AND  SLOPES  OF  THE  NODAL  POINTS 

PHI  =  STATE  TRANSITION  MATRICIES  FOR  EACH  MODE 

GAMMA  =  INPUT  TRANSITION  MATRIX 

A  =  DIAGONAL  MATRIX  CONSISTING  OF  PHI 

B  =  INPUT  MATRIX  OF  GAMMA  AND  CONTROL  NODE  SLOPES 

BN  =  NOISE  INPUT  MATRIX  OF  GAMMA  AND  NOISE  NODE  SLOPES 

DAMP  =  DAMPING  FACTOR 

SAMPT  =  SAMPLING  TIME 

IMPLSE  =  IMPULSE  INPUT  FUNCTION 

TCX,  TCY,  TCZ  =  CONTROL  TORQUE  VALUES 

IMPX,  IMPY,  IMPZ  =  AXIS  IMPULSE  NOISE  VALUES 

ENERGY  =  SYSTEM  ENERGY  COST  VALUE  FOR  A  GIVEN  POINT  IN  TIME 

CNTCST  =  SYSTEM  CONTROL  COST  VALUE  FOR  A  GIVEN  POINT  IN  TIME 

COST  =  TOTAL  SYSTEM  COST  VALUE  FOR  A  GIVEN  POINT  IN  TIME 

TOTCST  =  SYSTEM  COST  SUMMED  OVER  ALL  TIME 

MIN  =  NUMBER  OF  MINUTES  SYSTEM  WILL  BE  OBSERVED 

***Vf******  SAMPLE  OF  SPACEN  EXEC  FILE  ************** 

THIS  FILE  MUST  BEGIN  IN  COLUMN  1  AND  RUN  WITH  THE  FOLLOWING 
SEQUENCE  FOR  THE  INITIAL  RUN  OF  THE  PROGRAM: 

FORTVS2  SPACEN  (COMPILES  PROGRAM) 

SPACEN  (LOADS  AND  RUNS  PROGRAM) 

SUBSEQUENT  PROGRAM  RUNS  CAN  ELIMINATE  "FORTVS2  SPACEN"  IF  NO 
CHANGES  HAVE  BEEN  MADE  TO  THE  PROGRAM. 

CP  DEF  STOR  2M 

FI  4  DISK  KLAMA  OUTPUT  B  (PERM 
FI  30  DISK  XI  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 
FI  31  DISK  MODENG  SPACEN  A  (RECFM  F  BLOCK  80  PERM 

FI  32  DISK  TORQUE  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 

FI  33  DISK  ENERGY  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 

FI  34  DISK  MDECST  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 

FI  35  DISK  COUNT  INPUT  A  (RECFM  F  BLOCK  80  PERM 
FI  40  DISK  UTILITY  DATA  A  (RECFM  F  BLOCK  80  PER 
FI  41  DISK  RUN  DATA  A  (RECFM  F  BLOCK  80  PERM 

FI  42  DISK  KUG69  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 

FI  43  DISK  KUG23  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 

FI  44  DISK  KUG55  OUTPUT  A  (RECFM  F  BLOCK  80  PERM 

LOAD  SPACEN 
START  *  NOXUFLOW 


*********************************************************** 

PI  =  4. ODO  *  ATAN(1. ODO) 

SAMPT  =  0.  0 
DAMP  =  0.  0 
MODAL  =  0 
IMPX  =  0.  ODO 
IMPY  =  0.  ODO 
IMPZ  =  0. ODO 
NF  =  100 
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NG  =  100 
NH  =  100 
NZ  =  200 

Vr  V"  *  Vf  Vr  Vr  iV 


*****  number  of  minutes  the  system  will  be  observed  ****vV 

MIN  =  120.  0 

*****************iV*** ************************************** 

*****  SET  LENGTH  OF  MODAL  MODEL  ***** 

*■>>**1’?*******************************^'?********************** 


MODAL  =56 

*********************************************************** 
*****  read  lama  matrix  ***** 


READ(4,1001)  NAM 

READ( 4 , 1 002 ) ( LAMA( I ) , 1=1 , 100) 


*****  SCREEN  INTERACTION  ***** 

*********************************************************** 


V'f******iV**-.v*  STARTING  MODE  NUMBER  ************** 

SMODE  =  7 


V?  Vf  '.V  Vr  Vf  *  *  *  Vf  *  *  * 


NUMBER  OF  MODES  TO  SCAN  ************* 


MODE  =5 

EMODE  =  SMODE  +  MODE  -  1 

************  NOISE  INPUT  POSITION  ***'!’c*'>v*******yf*** 

NODE  =  55 
AXIS  =  1 

ynv***********  R  MATRIX  VALUE  **************** 


RM  =  lE-12 

********“.01  FOR  FULL  SAMPLING  TIME  ".05  FOR  REDUCED  ****** 


SAMPT  =  0. 01 
C 

Q  ***************  DAMPING  FACTOR  ************** 

DA.MP  =  0.  OOIDOO 
C 

DO  75  I  =  1,100 

READC42,1040)  ( UG69( I ,K) ,K=1 , 3) 

READ(43,1040)  ( UG23( I ,K) ,K=1 , 3) 
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READ(44,1040)  (UG55( I ,K) ,K=1 . 3) 
CONTINUE 

^-iV-.VVrVnViV********  BEGIN  RUN 


*************** 


DO  505  M  =  1,3 

WRITE  (41,700)  SMODE 
WRITE  (41,701)  MODE 
WRITE  (41,706)  EMODE 
WRITE  (41,702)  NODE 
WRITE  (41,703)  RM 
WRITE  (41,704)  SAMPT 
WRITE  (41,705)  DAMP 
WRITE  (41,707)  MIN 
WRITE  (41,708)  MODAL 

********  NOISE  AXIS  INPUT  AND  LOCATION 

IF(AXIS.  EQ.  DTHEN 
IMPX  =  1. ODO/SAMPT 
ELSEIF(AXIS.  EQ.  2)THEN 
IMPY  =  1. ODO/SAMPT 
ELSEIF(AXIS.  EQ.  3)THEN 
IMPZ  =  1. ODO/SAMPT 
ELSEIF(AXIS.  EQ.  4)THEN 
IMPX  =  1. ODO/SAMPT 
IMPY  =  1. ODO/SAMPT 
IMPZ  =  1. ODO/SAMPT 
END  IF 

COUNT  =  0 


********** 


*****VfrV*Vf** 


INITIALIZE  MATRICIES 


**************** 


DO  40  I  =  1,188 
DO  45  J  =  1,188 
PHI(I,J)  =  0. 0 
CONTINUE 
CONTI.NUE 

DO  60  I  =  1,188 
DO  65  J  =  1,3 
B(I,J)  =  0.0 
BN(I,J)  =  0. 0 
CONTINUE 
CONTINUE 

DO  70  K  =  7,100 
X1(K)  =  0.  0 
X2(K)  =  0. 0 
MODEN(K)  =0.0 
R.MODEN(K)  =0.0 
CONTINUE 


**'>V************')V***Vr*******yfVf**************,v**************,l 

***•-■*  BEGIN  MAIN  PROGRAM  ****,1 


45 


ODD 


*****  ESTABLISH  A,  B  AND  B"N0ISE"  MATRICIES  ***** 


DO  600  I  =  SMODE, MODAL 

WN  =  DBLE(SQRT(LAMA(I))) 

GMA  =  DAMP*WN/2. 0 
EGT  =  DEXP( -GMA*SAMPT) 

W1  =  DSQRT((WN**2)-(GMA**2)) 

COSWIT  =  DC0S(W1*SAMPT) 

SINT^’II  =  DSIN(W1*SAMPT) 

C 

IF(WN.  EQ.  0)THEN 

PHII(1,1,I)  =  EGT*C0SW1T 
PKII( 1,2,1)  =  SAMPT 
PHII(2,1,I)  =  0 
PHI  I  (2, 2, 1)  =  EGT’VCOSWIT 
C 

GAMMA(1,I)  =  0 
GAMMA(2,I)  =  0 

ELSE 

C 

PHII(1,1,I)  =  EGT’KCOSWIT  +  (GMA*(W1**( -1) )  )*SINW1T) 
PHII(1,2,I)  =  (W1**(-1))*EGT’'^SINW1T 
PHII(2,1,I)  =  -(WN**2)*(W1**(-1))*EGT’'^SIN’V1T 
PHII(2,2,I)  =  EGT->(C0SW1T  -  (GMA*(W1**( -1) ) )*SINW1T) 

C 

GAMMA(1,I)  =  (WN**(-2))*(1.  ODO-EGT*(COSW1T+(GMA/W1)*SINV1T) ) 
GAMMA(2,I)  =  (W1**(-1))*EGT^^SINW1T 
C 

ENDIF 

C 

600  CONTINUE 
C 

V  =  1 

c 

DO  610  K  =  SMODE, MODAL 
C 

PHI(V,V)  =  PHII(1,1,K) 

PHI(V,V+1)  =  PHII(1,2,K) 

PHI(V+1,V)  =  PHII(2,1,K) 

PHI(V+1,V+1)  =  PHII(2,2,K) 

C 

B(V,1)  =  GAMMA(1,K)*DBLE(UG69(K,1)) 

B(V,2)  =  GAMMA(1,K)*DBLE(UG69(K,2)) 

B(V,3)  =  GAMMA( 1,K)*DBLE(UG69(K,3)) 

B(V+1,1)  =  GAMMA(2,K)*DBLE(UG69(K,1)) 

B(V+1,2)  =  GAMMA(2,K)*DBLE(UG69(K,2)) 

B(V+1,3)  =  GAMMA(2,K)*DBLE(UG69(K,3)) 

C 

V  =  V+2 
C 

610  CONTINUE 
C 

DO  605  I  =  1,100 
UGVEX(I,1)  =0.0 
UGVEX(I,2)  =  0. 0 
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UGVEX(I,3)  =  0.  0 
IFCNODE.EQ.  23)  THEN 

UGVEXd.l)  =  UG23(I,1) 

UGVEX(I,2)  =  UG23(I,2) 

UGVEX(I,3)  =  UG23(I,3) 
ELSEIF(NODE.  EQ.  55)  THEN 
UGVEX(I,1)  =  UG55(I,1) 

UGVEX(I,2)  =  UG55(I,2) 

UGVEX(I,3)  =  UG55(I,3) 
ELSEIF(NODE.  EQ.  69)  THEN 
UGVEXdd)  =  UG69d,l) 

UGVEXdd)  =  UG69d,2) 

UGVEXdd)  =  UG69d,3) 

ENDIF 

605  CONTINUE 

C 

V  =  1 
C 

DO  620  K  =  SMODE, MODAL 
C 

BN(V,1)  =  GAMMAd,K)*DBLE(UGVEX(K,l)) 
BN(V,2)  =  GAMMAd,K)*DBLE(UGVEX(K,2)) 
BN(V,3)  =  GAMMAd,K)*DBLE(UGVEX(K,3)) 
BN(V+1,1)  =  GAMMA(2,K)*DBLE(UGVEX(Kd)) 
BN(V+1,2)  =  GAMMA(2,K)*DBLE(UGVEX(K,2)) 
BN(V+1,3)  =  GAMMA(2,K)*DBLE(UGVEX(K,3)) 

C 

V  =  V+2 
C 

620  CONTINUE 
C 

550  CONTINUE 
C 

c  ***********  establish  h,  f  and  r  matricies 

c 

DO  50  I  =  1,NH 
DO  55  J  =  1,NH 
Hd,J)  =  0.0 
Fd,J)  =  0.  0 
Gd,J)  =  0.0 
IFCI. LE. 3)THEN 

Ld,J)  =  0.0 

ENDIF 

55  CONTINUE 

50  CONTINUE 

C 

DO  61  I  =  1,3 
DO  66  J  =  1,3 

Rd,J)  =  0.  0 
RINVd,J)  =  0.0 
66  CONTINUE 

61  CONTINUE 
KQ  =  1 

DO  80  K  =  SMODE, EMODE 

H(KQ,KQ)  =  DBLE(LAMA(K)) 

H(KQ+1,KQ+1)  =  l.ODO 


********** 
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KQ  =  KQ+2 
80  CONTINUE 
C 

K  =  0 

DO  85  K  =  1,3 
R(K,K)  =  RM 
RINV(K,K)  =  l.ODO/RM 
85  CONTINUE 

C 

DO  88  I  =  1,2*M0DE 
DO  89  J  =  1,2*M0DE 
F(I,J)  =  PHI(I,J) 

89  CONTINTJE 

88  CONTINUE 

C 

C  *****  COMPUTE  G  MATRIX  AS  NEEDED  BY  RICDSD  SUBR  ***** 

C 

CALL  MATRAN(B,188,2*M0DE,3,BT,3) 

C 

CALL  MATMULC  B , 1 88 , 2*M0DE , 3 , RINV , 3 , 3 , TEMP , NH ) 

C 

CALL  MATMULC  TEMP , NH , 2*M0DE , 3 , BT , 3 , 2*M0DE , G , NG ) 

C 

C  BEGIN  RICCATI  GAIN  CALCULATIONS 

c 

CALL  RI CDSDC NF , NG , NH , NZ , 2*M0DE , 4*M0DE ,F,G,H,Z,W,ER,EI, WORK , 
+  SCALE,ITYPE,IPVS) 

C 

WRITE(41,*)  '  ' 

WRITE  (41,130) 

130  FORMAT  (/'  THE  CLOSED  LOOP  EIGENVALUES  ARE:  ’/) 

DO  140  I  =  l,2*MODE 

WRITE  (41,*)  ER(I),EI(I) 

140  CONTINUE 

WRITE  (41,150)  WORK(l) 

150  FORMAT  (/'  CONDITION  ESTIMATE  IS:  ' ,D26. 18) 

C 

C  ************  COMPUTE  GAIN  MATRIX  L  *********Vr*** 

C 

CALL  MATMULC  BT , 3 , 3 , 2*MODE , H , NH , 2*MODE , TEMP 1,3) 

C 

CALL  MATMUL(TEMPl,3,3,2*MODE,B,188,3,RR,3) 

C 

DO  103  I  =  1,3 

DO  104  J=l,3 

RR(I,J)=R(I,J)+RR(I,J) 

104  CONTINUE 

103  CONTINUE 

C 

CALL  DLINDS(3,RR,3,RRINV,3) 

C 

CALL  MATMULC RRINV , 3 , 3 , 3 , TEMPI , 3 , 2*MODE , BT , 3 ) 

C 
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noooo  Ci  n  ncionno  vononon 


CALL  MATMTJLC  BT ,  3 , 3 , 2*M0DE ,  PHI ,  188 , 2*M0DE ,  L ,  3 ) 

C 

WRITE(41,*)  '  ' 

WRITE(41,*)  '  GAIN  MATRIX  L  ' 

WRITE(41,*)  '  ' 

WRITE(41,*)  '  ROW  1  ROW  2  ROW  3' 

DO  9155  I=1,2*M0DE 

WRITE(41,1040)  (L(J,I),J=1,3) 

9155  CONTINUE 

WRITE(41,*)  ’  ’ 

*****  COMPUTATION  OF  TORQUES  ANTI  COSTS  ***** 


COUNT  =  0 
TOTCST  =  0.  ODO 
TIME  =  0.0 


*****************Vf****************************************Vf**** 

*****  SETS  LOOP  FOR  THE  NUMBER  OF  ITERATIONS  NECESSARY  ***** 
*****  TO  OBSERVE  THE  SYSTEM  FOR  DESIRED  LENGTH  OF  TIME  ***** 

*******5'.-***************Vc****Vf*****  ***************************** 


LOOP  =  INT((MIN*60. 0)/SAMPT) 

PRNT  =  INT(((MIN*60.  0)/SAMPT)/300.  0) 
PRNTG  =  INT(((MIN*60.  0)/SAMPT)/1000.  0) 

DO  200  N  =0,  LOOP 
TIME  =  DBLE(N)*SAMPT 


IF(N. EQ.  0)THEN 
IMPLSX  =  IMPX 
IMPLSY  =  IMPY 
I.MPLSZ  =  IMPZ 
ELSE 

IMPLSX  =  0. ODO 
IMPLSY  =  0.  ODO 
IMPLSZ  =  O.ODO 
END  IF 


**********  CONTROL  TORQUE  EQUATIONS  *********** 

*******************y?********************yf****************** 


SUMl  =  O.ODO 
SUM2  =  0. ODO 
SUM3  =  0. ODO 
C 

DO  210  CT  =  1,  MODE 

CTADJ  =  CT  +  (SMODE  -  1) 

SUMl  =  SUMl  +  L(1,2*CT-1)*X1( CTADJ)  +  L( 1 ,2*CT)*X2 (CTADJ) 

SUM2  =  SUM2  +  L(2,2*CT-1)*X1(CTADJ)  +  L(2,2*CT)*X2(CTADJ) 

SUMS  =  SUMS  +  L(3,2*CT-1)*X1(CTADJ)  +  L( 3 ,2*CT)*X2(CTADJ) 

210  CONTINUE 
C 
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ooooo  oonoonoorjoooooooooooooonnnooooooooooon  onooooo 


************ 


EVALUATION  AT  ZERO  CONTROL 


*********Vf** 


IF(M.  EQ.  DTHEN 
TCX  =  0.  ODO 
TCY  =  0. ODO 
TCZ  =  0.  ODO 
ELSE 

TCX  =  SUM1*(-1.  ODO) 

TCY  =  SUM2*(-1.  ODO) 

TCZ  =  SUM3*(“1. ODO) 

ENDIF 

IF(N.  EQ.  0)THEN 

WRITE  (32,*)  'IMPULSE  X  AXIS,  IMPULSE  Y  AXIS,  IMPULSE  Z  AXIS' 
WRITE  (32,*)  '  ' 

WRITE  (32,1040)  IMPLSX,  IMPLSY,  IMPLSZ 
WRITE  (32,*)  '  ' 

WRITE  (32,*)  'CONTROL  TORQUES  TCX,  TCY,  TCZ' 

WRITE  (32,*)  '  ' 

ENDIF 

IF  (N. LE.  50)  THEN 

WRITE  (32,2000)  TIME,  TCX,  TCY,  TCZ 
ENDIF 

IF  (MOD(N,PRNTG). EQ. 0)  THEN 

WRITE  (32,2000)  TIME,  TCX,  TCY,  TCZ 
ENDIF 

IF  (N. LE.20)  THEN 

WRITE(30,1035)  TIME ,X1( 1) ,Xl( 5 ) ,X1( 10) ,X1( 15 ) ,X1( 20) 
WRITE(31,1035)  TIME ,X2( 1) ,X2C 5) ,X2( 10) ,X2( 15) ,X2( 20) 

ENDIF 

IF  (M. EQ. 1)  THEN 
IF  (MOD(N,PRNTG).EQ.  0)  THEN 

WRITE(30,1036)  TIME ,X1(9) ,X1( 10) ,X1( 11) ,X1( 12) ,X1( 13) ,X1(14) 
WRITE(30,1036)  TIME,X1( 1) ,X1(4) ,X1( 14) ,X1(24) ,X1( 34) ,X1(44) 
WRITE(31, 1035)  TIME,X2( 1) ,X2(5) ,X2( 10) ,X2( 15) ,X2(20) 

ENDIF 

ELSE 

IF  (MOD(N,PRNT).EQ. 0)  THEN 

WRITE(30,1036)  TIME ,X1( 1) ,X1(4) ,X1( 24) ,X1(44) ,X1( 74) ,X1( 94) 
WRITE! 30, 1036)  TIME,X1( 1) ,X1(4) ,X1( 14) ,X1(24) ,X1(34) ,X1(44) 
WRITE! 31, 1035)  TIME,X2( 1) ,X2(5) ,X2( 10) ,X2( 15) ,X2(20) 

ENDIF 

ENDIF 


IF  (MOD(N,PRNTG).EQ.  0)  THEN 
COUNT  =  COUNT+1 
ENDIF 

*********************************************************** 
********  SYSTEM  COST  FUNCTION  CALCULATION  ******** 

**iV-VVf  ************■>’?*************■**************************** 
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nnnno  oooonoooo 


SUMC  =  0.  ODO 
ENERGY  =  0. ODO 
CNTCST  =  0.  ODO 
COST  =  0.  ODO 
C 

DO  230  CF  =  7, MODAL 

MODEN(CF)  =  MODEN(CF)+(Xl(CF)**2)*LAMA(CF)+X2(CF)**2 
SUMC  =  SUMC+(X1(CF)**2)*LAMA(CF)+X2(CF)**2 
230  CONTINUE 

ENERGY  =  SUMC 

CNTCST  =  (TCX**2)*RM+(TCY**2)*RM+(TCZ**2)*RM 
COST  =  ENERGY  +  CNTCST 
TOTCST  =  TOTCST  +  COST 

IF  (MOD(N,PRNTG).EQ.  0)  THEN 

WRITE(33,2000)  TIME .ENERGY, CNTCST, COST 
ENDIF 

IF(N.  GE.  (LOOP-50))THEN 

WRITE  (34,3002)  MODE, TOTCST 
ENDIF 

*****vr**  STATE  UPDATE  EQUATIONS  ********** 

*****'3VVr***Vf*******Vr***************Vc***************‘>V******** 

DO  220  KA  =  7, MODAL 
K  =  KA-6 
C 

X1T=PHII(1,1,KA)*X1(KA)+PHII(1,2,KA)*X2(KA)+B((2*K-1) ,1)*TCX+ 
+  B( ( 2*K-1) ,2)*TCY+B( ( 2*K-1) ,3)*TCZ+BN( ( 2*K-1) , 1)*IMPLSX+ 

+  BN( ( 2*K- 1) , 2)*IMPLSY+BN( ( 2*K- I ) , 3)*IMPLSZ 

C 

X2T=PHI I ( 2 , 1 , KA ) *X 1 ( KA ) +PHI I ( 2 , 2 , KA ) *X2 ( KA ) +B ( 2*K , 1 ) *TCX+ 

+  B ( 2*K , 2 )*TCY+B( 2*K , 3)*TCZ+BN( 2*K , 1 )*IMPLSX+BN( 2*K , 2 )* 

+  IMPLSY+BN(2*K,3)*IMPLSZ 

C 

Xl(KA)  =  XIT 
X2(KA)  =  X2T 
C 

220  CONTINUE 

C 

200  CONTINUE 
C 

WRITE  (35,3000)  COUNT 
RTOTAL  =  TOTCST 

C  WRITE  (34,3002)  MODE, RTOTAL 
WRITE  (34,3002)  MODE, TOTCST 
C 

DO  235  K  =  7, MODAL 

RMODEN(K)  =  MODEN(K) 

WRITE  (31,3002)  K,  RMODEN(K) 

235  CONTINUE 
C 

666  CONTINUE 
C 
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* 

IFCM.  EQ.  DTHEN 

* 

M0DE=10 

* 

EMODE=16 

* 

ELSEIFCM.  EQ.  2)THEN 

* 

MODE=15 

* 

EMODE=21 

* 

ELSEIFCM.  EQ.  2)THEN 

* 

MODE=20 

* 

EMODE=26 

* 

ELSEIFCM.  EQ.  2)THEN 

* 

MODE=25 

* 

EMODE=31 

* 

ELSEIFCM.  EQ.  3)THEN 

* 

MODE=30 

* 

EMODE=36 

* 

ELSEIFCM.  EQ.  4)THEN 

* 

MODE=35 

* 

EM0DE=41 

* 

IFCM.  EQ.  DTHEN 

* 

M0DE=40 

* 

EM0DE=46 

* 

ELSEIFCM.  EQ.  2)THEN 

* 

MODE=45 

* 

EMODE=5 1 

* 

ELSEIFCM.  EQ.  3)THEN 

* 

MODE=50 

* 

EMODE=56 

* 

END  IF 

*505 

CONTINUE 

c 

c 

c 

c 

c 

700 

701 

702 

703 

704 

705 

706 

707 

708 
1001 
1002 

1004 

1005 
1008 
1010 

1035 

1036 
1040 
1050 
2000 
2001 
3000 


*****  FORMAT  STATEMENTS  ***** 

*?’?****************************■)’«>**************************■>> 


(' 

(' 

(' 

(’ 

(’ 

(’ 

(' 

(' 

(' 


FORMAT 
FORMAT 
FORMAT 
FORMAT 
FORMAT 
FORMAT 
FORMAT 
FORMAT 
FORMAT 
FORMATC IX, 
FORMATdX, 
FORMATC IX, 
FORMATC IX, 
FORMATC IX, 
FORMATC Al) 
FORMATC'  ' 
FORMATC'  ' 
FORMATC'  ' 
FORMATC'  ' 
FORMATC IX. 
FORMATC' 
FORMATC 14) 


,' STARTING  MODE  NUMBER;  ’,12) 

,' NUMBER  OF  MODES  SCANNED;  ',12) 
input  NODE;  ',13) 

,' INITIAL  R  VALUE;  ',E12.4) 

,' SAMPLING  TIME;  ' ,E12. 4) 

, 'DAMPING  FACTOR;  ' ,E12. 4) 

,'LAST  CONTROLLED  MODE;  ’,12) 

,' OBSERVATION  TIME;  ',F5.1,'  MINUTES') 
,'SIZE  OF  MODAL  MODEL;  ’,13,'  MODES') 
A6) 

8E15.  8) 

//) 

60X,E11. 5) 

////) 

,F7.  2,2X,5CE12. 6,2X)) 

,F7. 2,1X,6CE11. 5,1X)) 

,3CE15.  8,5X)) 

,4CE12.  6,2X)) 

F7. 2,3X,3CE15.  8,3X)) 

,T5,E15.  8) 
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3001  FORMAT(F7. 2,2X,E12. 5) 

3002  FORMAT(I3,2X,E12. 5) 

C 

599  STOP 
END 
C 

Q  Vci'r******************^'-************************************** 

C  SUBROUTINE  TO  MATRIX  MULTIPLY 

Q  tViViViV^oV***************************************************** 

C 

SUBROUTINE  MATMUL(M1,LD1,R1,C1,M2,ID2,C2,MP,LD3) 

INTEGER  R1,C1,C2,LD1,LD2,LD3 
REAL*8  M1(LD1,1'),M2(LD2,1),MP(LD3,1),SUM 
C 

DO  650  I  =  1,R1 

DO  660  J  =  1,C2 
SUM  =  0. ODO 
DO  670  K  =  1,C1 

SUM  =  SUM+M1(I,K)*M2(K,J) 

670  CONTINUE 

MP(I,J)  =  SUM 
660  CONTINUE 

650  CONTINUE 

RETURN 
END 
C 

C  SUBROUTINE  TO  TRANSPOSE  A  MATRIX 

Q  yf^VV<'^VVf>VVfVf^VV?VrVf*^VVfVfVfVc’V?*Vr7V*VfiV*Vr'>V‘5VVc^Vf'>VVf*iV‘>V'>V^V^V^V‘iV'>V'>VTV^V*>V’*Vf^V^V^V>Vyr>V':^r^V>V^V 

c 

SUBROUTINE  MATRAN( MX , LDX , R 1 , C 1 , MT , LDT) 

INTEGER  Rl, Cl, I, J, LDX, LDT 
REALMS  MX(LDX,1) ,MT(LDT,1) 

C 

DO  680  I  =  1,R1 

DO  690  J  =  1,C1 

MT(J,I)  =  MX(I,J) 

690  CONTINUE 

680  CONTINUE 

RETURN 
END 
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APPENDIX  B.  KARHUNEN-LOEVE  MODE  PROGRAM 


This  program  computes  the  Karhunen-Loeve  mode  shapes  as  discussed  in  Chapter 
IV.  Subroutine  GENICAA  is  called  to  generate  the  covariance  matrix  of  the  modal  co¬ 
ordinates.  RAA  =  ECa(t)  a‘''(t)]]  .  One  subroutine  from  the  IMSL  math  librar\' 
(DMRRRR)  is  used  for  matrix  multiplication  and  two  EISPACK  subroutines  are  used 
to  compute  the  eigenvalues  and  eigenvectors  of  the  covariance  matrix. 

The  subroutine  GEXRAA  as  noted  above  generates  the  covariance  matrix,  RAA, 
of  the  modal  coordinates  of  the  space  station  discussed  in  Chapter  I.  The  algorithm: 

•  Computes  the  natural  frequency,  damping  coefficient,  y,  and  damped  frequency  of 
each  mode 

•  Computes  the  variances  and  covariances  based  on  these  parameters  as  given  by 
equations  19,  20,  21 

The  KL  mode  shapes  are  found  by  solving  the  eigenvalue  problem,  equation  12. 


R.AA*Q  (f)  —  ?.<!) 


where  the  matrix  Q  is: 


Q 


=  diag 


2 

CO7  0 

0  1 


cUjj  0 
0  1 


CO56  0 
0  1 


(35) 


(36) 


and 

•  /  is  the  eigenvalue  of  the  product,  RAA  *  Q 

•  (f)  represents  the  orthonormal  eigenvectors 

R.-\.‘\  "  Q  is  not  a  symmetric  matrix  so  the  eigenvalue  problem  is  changed  to: 

[Q'  =  R.^AQ''q  [Q’ / [Q”0]  (37) 


Let 


1/  =  Q'  and  0  =  Q-'  2 


and  finally: 

[q"raaq"]i/  =  ;.i? 


(3S) 


(39) 
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The  EISPACK  routines  for  symmetric  matrices  are  used  to  solve  the  eigenvalue  prob¬ 
lem.  equation  39,  and  equation  38  is  used  to  compute  the  KL  mode  shapes. 

This  program  is  run  by  executing  the  following  EXEC  file: 

•  FORTVS2  KLMODES2 

•  FI  4  DISK  THESIS  INPUT  B  (PERM 

•  FI  6  DISK  RAA  OUTPUT  T  (PERM 

•  FI  8  DISK  EIGVEC  OUTPUT  T  (PERM 

•  FI  10  DISK  EIGVAU  OUTPUT  T  (PERM 

•  FI  09  DISK  EIGEN  ERROR  B  (PERM 

•  FI  12  DISK  NATMODES  DOOl  B  (PERM 

•  FI  13  DISK  KUAMAQ  OUTPUT  B  (PERM 

•  FI  14  DISK  R.AA  OUTPUT  (PERM 

•  FI  17  DISK  NERROOl  DATA  A  (PERM 

•  P  DEE  STOR  2M 

•  EXEC  M.ATHP.ACK  EISPACK 

•  EXEC  TDISK  5  DIS 

•  UOAD  KLMODES2  (START 


PROGRAM  KLMODE 

VfVf-.'fi'n'fiVVf.v****  Var iabl®  definitions***************************’*^*****’''' 
lama  =  vector  of  natural  frequencies  * 

ugvex=  matrix  of  the  modal  amplitudes  and  modal  slopes  * 

natraod=vector  containing  the  natural  mode  shapes  that  correspond  * 
to  the  KL  mode  shapes  * 

*  klama=vector  of  natural  frequencies  ordered  according  to  th  KL  * 

*  mode  shapes  * 

*  klmod=vector  containing  the  largest  eigenvectors  in  column  of  the  * 

eigenvector  matrix  eigenvec  * 


*****■>'? 

* 

"it 
it 


CHARACTER*6  NAM 

DOUBLE  PRECISION  LAMA( 100 , 1) ,UGVEX( 684 , 100) ,D( 100) ,E( 100) 

DOUBLE  PRECISION  KLM0D( 100) ,Y,N1 ,Q( 100 , 100) 

DOUBLE  PRECISION  GAMA, DAMP, W(50) ,L( 50) ,M(50) ,B2( 1,50) ,Q2(100,100) 
DOUBLE  PRECISION  K( 1 ,50) ,A,B ,C ,C0EF,RAA( 100 , 100) ,APB ,APC ,TEMP1 
DOUBLE  PRECISION  TEMP2 ,TEMP3,TEMP4 ,TEMP5 ,TEMP6 ,Q2RAA( 100 , 100) 
DOUBLE  PRECISION  EIGVEC( 100, 100) ,Z( 100, 100) ,RAAQ( 100, 100) 

DOUBLE  PRECISION  Q2RQ2( 100 , 100) ,Q2INV( 100 , 100) ,KLAMA( 50) 

INTEGER  KK,LL,ZZ,MODR( 100) ,M0DC( 100) ,MM,NATM0D( 100) ,NUM 

DO  160  1=1,100 

DO  170  J=l,100 

Q2(I,J)=0. ODOO 

Q2INV(I,J)=0.  ODOO 

RAA(I,J)=0. ODOO 
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Q(I,J)=0.  ODOO 
170  CONTINUE 
160  CONTINUE 

READ(4,1001)  NAM 

READ(4,1002)  (LAMA( I , 1) , 1=1 , 100) 

READ(4,1001)  NAM 

READ(4,1002)  ( (UGVEX( I , J) , 1=1 ,684) , J=1 , 100) 

C  NODE  IS  THE  DISTURBANCE  INPUT  LOCATION  . 

C  AND  DILOC  IS  THE  COMPUTED  LOCATION  OF  THE  MODAL  SLOPE  IN  UGVEX 
NODE=69 

DISLOC  =  4+6*(N0DE-l) 

C  NUMODE  IS  THE  NUMBER  OF  THE  LAST  MODE  CONSIDERED 
NUMODE  =56 

CALL  GENRAA(NUMODE, DISLOC, LAMA, UGVEX, W,L,M,B2,K,RAA) 
WRITE(16,1005)W(1),L(1),M(1),B2(1,1),K(1,1) 

LL=2 

KK=1 

DO  150  J=7, NUMODE 
Q2(KK,KK)=DSQRT(LAMA(J,1)) 

Q2(LL,LL)=1.  ODOO 
Q2INV(KK,KK)=1/Q2(KK,KK) 

Q2INV(LL,LL)=1. ODOO 

KK=KK+2 

LL=LL+2 

150  CONTINUE 

CALL  DMRRRRC 100 , 100 ,Q2 , 100 , 100 , 100,Q2 , 100 , 100 , 100 ,Q, 100) 

CALL  DMRRRRC 100 , 100 ,RAA, 100, 100 , 100 ,Q, 100, 100, 100, RAAQ, 100) 

DO  15  1=1,100 
WRITE(14,1007)  RAA(I,I) 

WRITE(15,1007)  Q2(I,I) 

15  CONTINUE 

CALL  DMRRRRC 100 , 100 , Q2 , 100 , 100 , 100 , RAA , 100 , 100 , 100 , Q2RAA , 100 ) 
CALL  DMRRRRC 100 , 100 , Q2RAA , 100 , 100 , 100 ,Q2 , 100 , 100 , 100 ,Q2RQ2 , 100) 
CALL  TRED2C100,100,Q2RQ2,D,E,Z) 

CALL  TQL2C100,100,D,E,Z,IERR) 

IFCIERR.  NE. 0)THEN 
WRITEC9,*)'ERROR=' ,16 
C  PAUSE 

END  IF 

CALL  DMRRRRC 100 , 100 ,Q2INV, 100 , 100 , 100 , Z , 100 , 100 , 100 ,EIGVEC , 100) 
WRITEC6,1005)CCRAACI,J),  J=1 , 100) , 1=1 , 100) 

WRITEC 8, 1006) CCEIGVECC I, J),J=1, 100) ,1=1,100) 

WRITEC 10, 1007)CDC I) ,1=1,100) 

MM=1 

DO  10  1=100,1,-1 
KLM0DCMM)=0. ODOO 
DO  20  J=l,100 

IFCDABSCKLMODCMM)).  LT.  DABSCt.IGVECC  J, I ) ) )THEN 
KLMODCMM)=EIGVECCJ,I) 

MODRCMM)=I 

MODCCMM)=J 

Y=CDBLECM0DCCMM)))/2. ODOO 
I1=INTCY) 

N1=Y-CDBLECI1)) 

IFCDABSCNl).  NE.  0.  ODOO)THEN 
NATM0DCMM)=CM0DCCMM)+l)/2  +  6 
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ELSE 

NATM0D(MM)=M0DC(MM)/2  +  6 

ENDIF 

END  IF 

20  continlt; 

MM=MM+1 

10  CONTINUE 

DO  30  1=1,100 

WRITE ( 12,100)  MODR( I ) ,MODC( I) ,NATMOD( I) ,KLMODC I) 

30  CONTINUE 

ZZ=50 

DO  180  1=1,100,2 

KLAMAC  ZZ )=LAMA( NATMODC I ) , 1) 

180 

WRITEC 13 , 1005) (KLAMAC I) , 1=1,50) 

NUM=7 

DO  11  1=1,99,2 

IF( DABS( KLMOD( I ) ) .  GT.  DABS(KLMOD( I+l) ) )THEN 
WRITE ( 1 7 , 200 )  NUM , DSQRT( 2*( 1 -DABS( KLMOD( I ) ) ) ) 

ELSE 

WR I TE ( 1 7 , 2 0 0 )  NUM , DSQRT( 2* ( 1 -DABS ( KLMOD ( I + 1 ) ) ) ) 

ENDIF 

NUM=NUM+1 

11  CONTINUE 

100  F0RMAT(1X,I3,2X,I3,2X,I3,2X,E15. 8) 

200  F0RMAT(1X,I2,2X,E15. 8) 

1001  F0RMAT(1X,A6) 

1002  F0RMAT(1X,8E15. 8) 

1005  F0RMAT(1X,5E15. 8//) 

1006  F0RMAT(1X,5E15. 8//) 

1007  F0RMAT(1X,5E15. 8//) 

STOP 

END 

SUBROUTINE  GENRAA  (NUMODE,DISLOC,LAMA,UGVEX,W,L,M,B2,K,RAA) 
DOUBLE  PRECISION  LAMA( 100 , 1) ,UGVEX( 684 , 100) 

DOUBLE  PRECISION  GAMA, DAMP, W( 50) ,L(50) ,M(50) ,B2( 1,50) 

DOUBLE  PRECISION  K( 1 , 50) , A,B ,C ,COEF,RAA( 100 , 100) , APB , APC ,TEMP1 

DOUBLE  PRECISION  TEMP2 ,TEMP3,TEMP4,TEMP5 ,TEMP6 

DAMP=0. OOIDOO 

GAMA=DAMP/2. ODOO 

DO  100  MODE=7,NUMODE 

W(MODE-6)=DSQRT(LAMA  (MODE,l)) 

L(MODE-6)=W(MODE-6)  '-DSQRT(  1.  0D00-GAMA*GAMA) 

M( MODE - 6 )=GAMA*W( MODE - 6 ) 

100  CONTINUE 

DO  110  J=7,NUMODE 

B2(l,J-6)=UGVEX(DISLOC+4*(J-l),J) 

K(l,J-6)=B2Cl,J-6)/L(J-6) 

C  K(l,J-6)=UGVEX(DISLOC,J)/L(J-6) 

no  CONTINUE 

DO  130  1=1,50 
DO  140  J=l,50 
A=M(I)+N(J) 

B=L(I)-L(J) 

C=L(I)+L(J) 
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C0EF=K(1,I)*K(1,J) 

APB=2*(A*A+B*B) 

APC=2*(A*A+C*C) 

RAA(  2*1 - 1 , 2*J- 1 )=COEF*( A/APB-A/APC ) 
TEMP1=L(J)*B/APB+L(J)*C/APC 
TEMP2=M( J)*A/APB-A*M( J)/APC 
RAA (2*1-1, 2* J ) =COEF* ( TEMP 1 -TEMP2 ) 

TEMP3=( L( I )*L( J)*A-L( I )*M( J)*B-L( I )*M( I )*B+M( I )*M( J)*A)/APB 
TEMP4=(L(I)*L(J)*A-L(J)*M(J)*C-L(I)*M(I)*C-M(I)*M(J)*A)/APC 
RAA( 2*1 , 2*J)=C0EF*(TEMP3+TEMP4) 

TEMP5=L( I )*C/APC-L( I )*B/APB 
TEMP6=M( I )*A/APB-M( I )*A/APC 
RAA(  2*1 ,2*J-1)=C0EF*(;TEMP5-TEMP6) 

140  CONTINUE 
130  CONTIN^'E 
RETURN 
END 
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